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o ■ 

' We perform a numerical test of a relativistic heavy quark(RHQ) action, recently proposed by 

Tsukuba group, in quenched lattice QCD at a ~ 0.1 fm. With the use of the improvement parameters 
previously determined at one-loop level for the RHQ action, we investigate a restoration of rotational 
symmetry for heavy-heavy and heavy-light meson systems around the charm quark mass. We focused 
on two quantities, the meson dispersion relation and the pseudo-scalar meson decay constants. It is 
OO . shown that the RHQ action significantly reduces the discretization errors due to the charm quark 

mass. We also calculate the S-state hyperfine splittings for the charmonium and charmed-strange 
' mesons and the Ds meson decay constant. The remaining discretization errors in the physical 

^ , quantities are discussed. 

m 
cn 

I. INTRODUCTION 

\^ • For a search of new physics beyond the standard model through flavor physics, a precise determination of physical 
' quantities such as quark masses and hadronic matrix elements associated with heavy mesons is required within a few 
% accuracy. Although, in principle, lattice QCD calculation is an ideal tool for this purpose, it suffers from large 
discretization effects due to the charm and bottom quark masses: ar7ic>0.3, amt,>l in quenched approximation and 
JL 1 amc>0-h, aTOb>1.5 in unquenched simulation with current computational resources. If we adopt the 0(a) improved 
Wilson quark action for the heavy quarks, the leading cutoff errors are expected to be 0((aTOg)"). In order to achieve 
(-H ' a few % accuracy, it is necessary to reduce the discretization errors for the heavy quarks to the same level for the 
• light quarks, which is ©((oAqcd)^)- 
.5^ _ For this end an on-shell ©(oAqcd) improved RHQ action has been proposed in Ref. 1], which extends the well 
^ . known on-shell improvement program 0, 0, S S 0] to massive quarks with amq ^ 0{1). The action works better as 
^ ' the lattice spacing becomes smaller. This is a fascinating feature from a view point of controlling the systematic errors 
5^ , coming from the cutoff effects. Another important point is that the action allows to treat the charm and bottom 
quarks simultaneously. In case of NRQCD, widely used for calculation of the bottom quark physics, on the other 
hand, it is theoretically impossible to take the continuum limit and difficult to treat the charm quark. 
The explicit form of the RHQ action is given by 



^RHQ ^ ^ 



moq{x)q{x) -|- q{x)j4:D4q{x) -(- iy'^q{x)jiDiq{x) 
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where four improvement parameters v, r^, ce and cb are relevant, while Vt is redundant. The leading cutoff effects 
of 0{{amQ)^) in this formulation can be removed by adjusting v and quark field renormalization factor as a function 
of aruQ and gauge couphng constant g. We can also remove the next-to-leading cutoff effects of 0{{amQ)"- a Aqcd) 
by adjusting r^, ce and cb- Note that it is recently pointed out that one can remove 0{{amQ)"' a Aqcb) errors in the 
spectral quantities such as masses by adjusting only two parameters, cb = ce and I'Q- However it is still true that 4 
parameters, v, rg, ce and cb, are necessary to remove all 0{{amQ)"' a Aqcd) errors in on-shell matrix elements. Once 
these are achieved, we are left with only the cutoff effects of 0{f2{amQ) {a Aqqd)^) where f2{a.mQ) is an analytic 
function of amq around amg = 0. We assume f2{amQ) ^ 0{1) for the massive quarks with mqa ^ 0(1). In the 
massless limit v and become unity and ce and cb agree with csw- Since, at present, it is difficult to determine Vs, 
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cb nonperturbatively, we have performed a perturbative determination of the four parameters at one- loop level in 
a mass dependent wayQ. In this case the remaining leading cutoff effects are 0{al fQ'^\amQ)) for the RHQ action. 
Similarly, we have also determined the mass-dependent renormalization constants and the improvement coefficients 
for the vector and axial- vector currents at one- loop level 9]. 

In this work we study the discretization effects of the RHQ action with the perturbatively determined improvement 
parameters using two different gauge actions in quenched approximation at a finite lattice spacing of a ~ O.lfm. We 
employ four heavy quark masses around the charm quark mass. In order to investigate the discretization effects, we 
calculate the dispersion relation both for the heavy-heavy and heavy-light mesons and the space-time symmetry of 
the pseudo-scalar meson decay matrix elements. For comparison, the same calculation is repeated with the heavy 
clover quark action. We observe sufficient reduction of (otoq)" errors in the RHQ action. In addition we extract the 
charmed-strange and charmonium hyperfine splittings and the Ds meson decay constant fjj^ at the physical charm 
quark mass. We compare our results with previous ones and discuss the cutoff effects on these quantities with the 
RHQ quark action. 

This paper is organized as follows. In Sect. [Til we explain the RHQ action and fix the notations. Simulation details 
are summarized in Sect. IIIII In Sect. IIVI we give a comparison between the RHQ action and the clover action by 
showing the dispersion relation and the space-time symmetry. We also present the results of physical quantities such 
as the hyperfine splittings and the decay constants, and discuss their cutoff errors in detail. A comparison of our 
results with the previous ones for the hyperfine splittings and the decay constants is also shown. In Sect. |V] we give 
our conclusion. A brief review of recent works related to this formulation can be found in Ref . [lOj . 



II. FORMULATION 



A. Actions 

For the gauge part we employ a renormalization-group (RG) improved gauge action proposed by Iwasaki[Tl| as well 
as the ordinary plaquette gauge action. For the quark part we use the clover quark actionQ for the light quarks and 
the RHQ action for the heavy quarks. We rewrite the RHQ action of Eq.([T]) to the following form with the use of 
hopping parameter k, which is more suitable for numerical simulations: 

^RHQ ^ j2qi^)D.M, (2) 



k=l,3 



-SxyCBu'^a^jFijix) ~ 5„jCEK^a4zF4i{x), (3) 

i<j i 

where the field strength in the clover terms is expressed as 

^)-Uh^)), (4) 

^X + l>,fl^XMJ (5) 

l>,fj.^x-fj.,iy^x-i^:t^' (6) 

Ul_-^_0^^Ux-fi-i>,fiUx-0,i/, (7) 

0,liUx+il-C',uUl^^. (8) 

As mentioned in the introduction, the parameters v, r^, cb and ce are already determined as a function of the heavy 
quark mass up to one-loop level. In this work we choose r* = 1. 
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B. Improvement of axial vector current 

The form of renormalized axial vector current with the ©(oAqcd) improvement is given in Ref.0]: 

A^ix) = y^q^/2^ZA,[q{x)-/^^5Qix)-C+^{q{x)A+j5Q{x)}~c^^{q{x^ 

+ci^{^»9(a;)}7<7M75Q(a;) - c2-^g(x)7^757,{AjQ(a;)} , (9) 

where q and Q denote the Ught and heavy quark fields, respectively. Za^^ is the finite renormalization factor connecting 
the lattice to the continuum MS scheme {Za = as defined in Ref. Q. Improvement coefficients for the 

temporal direction are in general different from those for the spatial direction: c^'f^'^ ^ '^^'^'^ ■ additional 
overall factor ^^k^^Ikq is associated with the field redefinition in the RHQ action of Eq.Q. Za^^ and c^^'^'^^ are 
calculated as a function of the quark masses amq and artiq up to one-loop level S] . With the aid of the equation of 

motion we can always choose — — 0. Covariant lattice derivatives A+ and A^ are defined as A+ = A^ + A^ 
and A^ — A ^ — A ^ , where A ^ and A ^ are lattice derivative acting on the right or left field as 

75 A^ Q{x) = q{x) 75 i[ U ^,{x)Q{x + fi) - Ul{x - fi)Q{x ~ fi) ], (10) 

g(x) A^ 75 (5(x) = ^[q{x + fi)Ul{x) - q{x - jl)U^{x - jl) ] '^5 Q{x). (11) 

Note that the definition of the improved current of Eq. ^ is slightly modified from that in Ref . . The difference will 
be explained in Appendix A. Hereafter we use a short-handed notation for the improvement terms such as 

0+ ^ -c+^{g-(a;)A+75Q(x)}, (12) 
O- EE -c^Jqix)A~j,Q{x)}, (13) 

= +c^JA,g-(a;)}7.7fe75Q(x), (14) 

0« = -c«j(x)7fe757.{A,Q(x)}. (15) 
The pseudo-scalar meson decay constant is defined as 

(OK(O)IP^(P)) = ipJpsiA^), (16) 

where \PS{p)) is the pseudo-scalar meson state with momentum p. The decay constant denoted as fpsiA^) identifies 
which component of the axial vector current is used. In the continuum limit /p5(v4^) for different ^ should agree 
with each other. 



III. SIMULATION DETAILS 
A. Simulation parameters 



We employ a single value of the gauge coupling constant, /? — 2.6 for the Iwasaki action and /3 = 6.0 for the 
plaquette action, on a L"^ x T = 24"^ x 48 lattice. Gauge configurations are generated by a 5-hit pseudo heat bath 
update supplemented by four over-relaxation steps. These configurations are then fixed to the Coulomb gauge at 
every 100 sweeps for the Iwasaki gauge action and 200 sweeps for the plaquette gauge action. We have accumulated 
300 gauge configurations for each gauge action with the RHQ action and 250 with the heavy clover quark action. 
With the use of the Sommer scale rg = 0.5fm the lattice cutoffs are determined as a^^(ro) = 2.0129(46) [GeV] [ll| and 
a^^('"o) — 2.1184(94)[GeV] [T^, respectively. The spatial lattice size in physical unit is approximately 2.4fm, which 
is large enough for the charmed mesons. 

Simulation parameters for the quark part are summarized in Table[T]for the Iwasaki and the plaquette gauge actions, 
where Mps/My represents the pseudo-scalar to vector meson mass ratio of the light-light and heavy-heavy mesons. 
For each gauge action we adopt three values of the light quark masses corresponding to Mps/My — 0.55 — 0.78 to 
cover the strange quark mass and four values of the heavy quark masses to sandwich the charm quark mass. 

For the light quarks we use the clover quark action with the nonperturbative value of the clover coefficient: Cg^ = 
1.50(5)pj|] for the Iwasaki action and 1.769p^] for the plaquette action. Here it is noted that Cg^ = 1.50(5) for the 
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Iwasaki action is taken from the preliminary result obtained in the infinite volume limit, which is 6 % larger than 
the final value Cg^ = 1.41 of Ref.p^] defined on a fixed physical volume. For the heavy quarks we adopt the RHQ 
action with the improvement parameters v, r^, ce and cb determined up to one-loop level with the mean-field(MF) 
improvement, details of which are explained in Appendix A. In order to remove O(aAQCD) errors at the massless 
point, we replace a massless part of c^; and cs by their nonperturbative value Cg^ as 

CE,B = c^'^s(aTOpoio) - c^'^B(aTOpoic = 0) + Cgw, (17) 

where the superscript PT represents the perturbative value up to one-loop level. 

For the light-light current we use nonperturbative values of the renormalization factor and the improvement coeffi- 
cients for the plaquette gauge action: = 0.807, hY = 1-28 and cY = 0.037[l3. For the Iwasaki gauge action we 
employ the mean-field improved values: Za = 0.86057, bA — 1.19998 and ca = 0.00864[l3. At present nonperturba- 
tive values are not available for this action. For the heavy-light and heavy-heavy currents, on the other hand, we use 

the mean-field improved values for c^^'^'^'^ and Za^ at the one-loop level (see Appendix A) for both gauge actions. 

In case of the plaquette gauge action we replace the massless part of c\ and Za , by the nonperturbative ones, 

and ZY- 

+ 2 -l-.PT/ \ 2 +.PT/ n^ I +,NP /1 o^ 

"^I,, ^ -9 "^A,, (a"^poic) - g (aTOpoio = 0) + c^' , (18) 

Za, - ZPJ(ampoie)-^5J(ampoic = 0) + Z^P. (19) 

In order to investigate a degree of improvement for the RHQ action we have made an additional simulation using 
the clover quark action both for the heavy and light quarks. Simulation parameters are given in Table [Til where we 
employ one value of the light quark mass and three values of the heavy quark masses roughly equal to lighter three 
for the RHQ action. 

On each configuration fixed with the Coulomb gauge, we invert the quark matrix employing the BiCGstab algorithm 
with the stopping condition that the residual must be smaller than 1.0 x 10"^*. For the heavy quarks we perform a 
fixed number of iterations. We choose 2T — 96 such that the stopping condition is always satisfied and it is assured 
that the heavy quarks can propagate from the origin to any point on the lattice. For both the light and heavy quark 
propagators we employ not only a local source but also an exponentially smeared source with a form of >lexp(— Br), 
where smearing parameters A and B are tuned to enhance an overlap with the ground state. Numerical values of A 
and B are listed in Table Hill for each combination of the gauge and quark actions. 



B. Measurement of two-point functions 

We measure the S-state (i.e. pseudo-scalar and vector) meson spectra for the light- light(L-L), heavy-light (H-L) and 
heavy-heavy(H-H) systems using the correlation functions projected onto zero spatial momentum state: 

X^(Oi(x,i)O^_j0)), (20) 

X 

where O — P ov V is understood. The subscripts S and L represent the smeared and local operators, respectively. We 
always adopt a local sink while taking both the local and smeared sources. Note that both the quark and anti-quark 
fields in Os are smeared. 

To extract the pseudo-scalar meson decay constant for the L-L, H-L, H-H systems, we calculate the correlation 
function 

Y,{^r\s,t)plm. (21) 

X 

where the superscript impr represents the 0{a) improved current given in Eq.(I5]). 

We also measure the meson correlation functions with finite spatial momenta given by 

ap= ^ X {(1,0,0), (1,1,0)}. (22) 

These correlation functions are used to calculate the dispersion relation of the S-state mesons and also to extract the 
decay constant using the temporal and spatial components of the axial vector current. 
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C. Fitting procedure 

The correlation functions in Eq. pO[) are expected to take the foUowing form for a large euclidean time separation: 
Ej(0t(f, i)O^(O)) = '^^^ e-aA/T/2 cosh(aM(T/2 - i)), (23) 



2aM 

7 . 7t 

I]^«P'^(f,t)pt(0)) = ^|__^e-«^^^/2sinh(aM(r/2-t)), (24) 



where M is a mass of the ground state allowed to couple to the operator. Matrix elements in the above expressions 
are given by 

Zps,. = (0|P5,l(0)|P5(p = 0)), (25) 
Zvs.. = (0|Fs,L(0)|y(p = 0)), (26) 
Z^.„,p. = (0KP'-(0)|P5(p = 0)), (27) 

where \PS(p = 0)) and \yij) = 0)) represent the pseudo-scalar and vector meson states at rest. We first extract 
aM by fitting the correlators of Eg. ([23|) . and then perform a fit of Eq. p4|) with aAI fixed. We employ the same 
fitting procedure for the correlation functions with finite spatial momenta. Since our statistics are not sufhcient to 
incorporate correlations between different time slices, we always use the uncorrelated fit for our analysis. We estimate 
statistical errors by the jackknife method with a bin size of 10 configurations to eliminate autocorrelations. In case 
that the correlated fit is possible, we use it to check the results obtained by the uncorrelated fit. We find that the 
results are consistent within statistical errors. 

The fitting ranges summarized in Table [TVl are chosen by investigating effective mass plots of the meson correlators 
presented in FiglU where we take K3 for the light quark and kq for the heavy quark as a representative case. Note that 
Kg roughly corresponds to the charm quark mass. We take similar fitting ranges for the correlators with finite spatial 
momenta, which are given in Tabl dlVI Figure [2] shows effective mass plots for the pseudo-scalar meson correlators 
with finite spatial momenta. 



IV. RESULTS 



A. Dispersion relation and space-time interchange symmetry 

In case that the improvement parameters are perturbatively determined up to one-loop level, the leading cutoff 

errors in the RHQ action is theoretically expected to be 0(aj/Q^'' (amg)), where f^\a'mQ) ^ 0(1) is assumed for 

amq ~ 0(1). We numerically check this theoretical expectation by investigating the dispersion relation of the S-state 

mesons and the space-time interchange symmetry for the pseudo-scalar meson decay constant. These quantities are 

(2) 

sensitive to the cutoff effects for the heavy quarks, and hence suitable to estimate a size of /q (otoq). 

We calculate an effective speed of light c^ff both for the pseudo-scalar and vector mesons by fitting the meson energy 
aE{p) as a function of the spatial momentum ap with the following form: 

{aE{p)f^clM' + ic^E{0)f- (28) 

In the continuum limit Cog should become unity. At finite lattice spacing, however, CcB deviates from unity due to 
the lattice cutoff errors. In Fig^{aE)^ is plotted as a function of (ap)'^, where the fitting result with Eq. is given 
by the solid line together with the continuum dispersion relation with CcS = 1 represented by the dashed line. We 
observe that the linearity of in is well satisfied and Ceff is close to unity. Fitted values of Ccff for the L-L, H-L 
and H-H cases are plotted in Fig. |4]for the pseudo-scalar mesons and in Fig. [5] for the vector mesons. Here it should 
be noted that in addition to finite quark mass errors CeS suffers from finite momentum corrections of 0{\ap\'^) so that 
Ceff could deviate from unity even for the massless quarks. Indeed Figd] shows that as the meson mass Ma decreases, 
Ceff becomes closer to unity within this uncertainty. In the heavy quark mass region around Ma ~ 1 — 2, Coff for the 
heavy clover quark action deviates from unity by about 7 — 10%. On the other hand, the RHQ action satisfies Coff = 1 
within 2 — 3% errors, which are comparable to the deviation for the L-L case. Since fitted values of Ccff for the vector 
mesons in Fig. [5] are consistent with those for the pseudo-scalar mesons within statistical errors, we use the values of 
Ccff determined from the pseudo-scalar meson dispersion relation in the following discussion. We observe no obvious 
difference in the results between the Iwasaki and plaquette gauge actions. 
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We also study the space-time symmetry of the pseudo-scalar meson decay matrix element defined by 

''-\Q\Af\PS{p)) \puV ^''^ 

where and Af' represent the spatial and temporal components of the renormalized axial vector current given in 
Eq.®. The pseudo-scalar meson state has finite spatial momentum of \ap\ = 27r/24. The ratio R is plotted in Fig[6] 
as a function of the meson energy Epg with the lowest finite spatial momentum for the L-L, H-L and H-H systems, 
where c^""""^^^ represents the partial replacement of the perturbative value for c\ by the nonperturbative one defined 

in Eg. p^ . while c^"^ means the perturbative value for c\ without this replacement. For the plaquette gauge action 
we employ Za^, defined in Eq.([T9l), though = 0.814 and = 0.807 agree with each other within 1%. Although 
because of the finite momentum corrections the ratio R could deviate from unity even for the massless quarks, it 
becomes consistent with unity within the statistical errors as the meson energy aE vanishes. For the massive mesons 
with aE ~ 1 — 2, on the other hand, the heavy clover quark action violates the space-time symmetry by about 7— 13%, 
while the RHQ action retains R = 1 within 6% errors. An intriguing observation is that the ratio R of the H-L system 
shows different aEps dependences between the Iwasaki and plaquette gauge actions: the ratio R decreases for the 
Iwasaki action as Eps increases, while it increases for the plaquette action. This different behaviors could come from 
a fact that the contributions of the 0(a) improvement operators are sizable for the plaquette action, whereas they are 
small for the Iwasaki action. This is observed in Figs 17] and [8] which show the relative contribution from each 0{a) 
improvement operator of Eas. (|12p - (jl5p to the axial-vector currents defined by 

j:,{A,{x,t)p{o)) ■ ^ ^ 

Dominant contributions always come from operators for the plaquette action, while their contributions are not so 
large for the Iwasaki action. In particular, this feature is more prominent for the H-L system. 

From the above analyses on Cog and R it can be concluded that the RHQ action succeeds in significantly reducing 
the {mQa)" errors in the heavy clover quark action. 



B. Physical quantities of S-state charmed mesons 

1. Physical points 

In order to obtain the meson spectra and the decay constants at the physical quark masses, we have to interpolate 
the heavy quark mass to the charm quark mass nic, while extrapolating the light quark mass to the w, d quark mass 
niud or interpolating it to the strange quark mass m^. Since we employ only 3 values of the light quark masses in 
our simulation, we consider only a linear extrapolation to the u, d quark mass. In the following the lattice spacing is 
always determined by the Sommer scale with vq = 0.5fm. 

The light-light pseudo-scalar meson masses are linearly fitted in 1/k as 

a'Mj,s = A 1-) , (31) 

where is determined from the vanishing point of [aMps]^ . K^d and Kg are determined so as to satisfy Mpg = 
= 135.0MeV and Mps = Mk = 497.7MeV, respectively. The fitting results of A and Kent are tabulated in 
Table El and Kud and Kg are given in Table IVTl 

We determine Kc in two different ways: matching Afpj''^ to Md^ — 1.9683GeV for the charmed-strange meson or 
My°^^ to Mji^ = 3.0969GeV for the charmonium, where the superscript pole represents a pole mass determined from 
an exponential fall-off of the meson correlator. Employing the following fitting functions 

aMlf = A + Bnuea^y + C«:La., + Dam';^^' (32) 

for the heavy-light meson masses with am^^''* — {l/nught — ^jKcrit)!'^ and 

aMP°'^ - A + BKheavy + C^La.y (33) 

for the heavy-heavy meson masses, we have determined two values of Kc-, which are given in Table IVTl 
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In order to estimate a magnitude of the cutoff errors, we also calculate the charmed meson spectra employing the 
kinetic mass defined by 

aM^''' = aMP°^ycls. (34) 

With the same fitting functions as Eqs. ip^ and (|33p we have also determined Kc{Ds, M"^'") and KciJ/ip, M^™) listed 
in Table IVTl From these results we observe that a difference of Kc between two physical inputs Mj/^ and Mjj^ is less 
than 0.2%, while a difference of Kc between MP°'° or M^™ is about 2%. In the following analysis, we always calculate 
all the physical quantities using both MP°'° and M^™, in order to estimate the systematic errors due to an ambiguity 
in the choice of MP°1'= or M^'". 



2. Hyperfine splitting for charmonium and charmed-strange meson 

Figure ini shows aMy dependence of the S-state charmonium hyperfine splitting aAM^ = aMy — aMpg, where 
X = pole or kin. In order to interpolate the results at the physical charm quark mass, we adopt the ansatz that the 
splitting is a polynomial of the inverse vector meson mass: 

aAM^ ^ A/{aM^) + B/{aM^f + CI{aM^f, (35) 

incorporating a property that the hyperfine splitting vanishes in the infinite quark mass limit due to the heavy quark 
symmetry. The interpolation lines are also plotted in FiglQl Using the fitting results for the parameters A, B and C 
given in Table|VlIl we obtain AM{J/^-t]c) in physical unit. AMP"^" at ndJ/i', MP°i<=) and AM'''" at Kc(J/'0, M^'") 
are tabulated in Table IVlTTl for each gauge action together with the experimental value. 

In FiglTO] we plot the S-state charmed-strange meson hyperfine splitting aAM^ = aMy — aMpg as a function 
of aMpg together with the interpolation lines which are obtained by employing the ansatz motivated by the heavy 
quark symmetry: 

aAM^ = {A + Bam^^s''')/{aM^g) + C/{aM^s?- (36) 

Using the fitting resuhs presented in Table HXl we obtain A M{D *-Ds) in physical unit. AMP"'^ at Hc{Ds, M^"^") and 
Ks, and am''™ at Kc(-Ds, Af'"") and hsted in Table IVUll for each gauge action together with the experimental 

value. 



3. Ds meson decay constants 

The heavy-light pseudo-scalar meson decay constant afps can be obtained from the temporal and spatial compo- 
nents of Eg. pB)) . In our calculation fps{A4) is determined from 

^aM^f afps = Z^.^p, / ^/ aMf,f (37) 

and fps{Ak) from 

a=^/2$fe = ^JaE^pf afps = Z aE^f / [lapk) . (38) 
where Z^impr and Z^impr are the decay matrix elements defined in Eq. (j27p . Note that only the improved axial vector 

4 k 

current with c^^^^^ is considered for the plaquette gauge action. In Fig[Tl]we plot a'^/^$p and a^/^$p as a function 
of l/(aMP°'°''''"). The interpolat ion lines are obtained by fitting the results with the following ansatz: 

a3/2$4^ ^ A + B/{aMlp''''') + C/{aMlf-'''y + Davi'^3ht^ (39) 
a3/2$^g = A + B/{aElf-''"^)+C/{aElf^'^y + Dam'^sht_ 

Using the fitted values of the parameters in Table IXllXIl we obtain /d, in physical unit. Table IXlTl lists the results of 
fps at Kc(I?s, AfP°'°) and Ks and fps at Kc{Ds,M^™) and Ks for each gauge action together with the experimental 
value. 
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C. Cutoff effects 

We now consider the cutoff effects in our results. Leading cutoff effects for the gauge part are 0(a^ Aq^q). The hght 
quark action also has 0{a^AQQj~,) errors, since the nonperturbative value of csw is employed for each gauge action. 

For the RHQ action, on the other hand, the leading cutoff effects are 0{as{^if' f^\a'mQ)) with as{^i) = [ji) / [At:) , 
which comes from the fact that the parameter v associated with the 0(1) kinetic term is only adjusted up to one-loop 
level. Since this error is responsible for the deviation of Cofr from unity, the mass dependence of Coff shown in FigsH] 
and [5] tells us that /q {niQa) is a smooth function of amg in the range of the heavy quark mass employed in our 

simulation. In addition, there exists the 0{as{n)'^g^\amQ)) errors originating from the heavy quark axial vector 
currents whose renormalization factors are determined up to one-loop level. These are the leading cutoff effects in the 
deviation of R from unity shown in Fig|Sl where we find fairly smooth amq dependence. 

Let us take into account these 0{as{^Y f^\amQ), as{fi)'^ g^p {airiQ)) effects in our error estimate using a difference 
of the charmed meson hyperfine splittings obtained with Mpoic and Mkin and also a difference of the charmed meson 
decay constants extracted from Af, and A4. For the hyperfine splittings we take the pole mass result as the central 
value and a difference between two results as a systematic error. In Table IXIIII our final result for the charmonium 
hyperfine splitting in physical unit is also presented, where the central value is AMpoic^J/tp — Vc), the first error is 
statistical and the second is a systematic error explained above. The second error, much larger than the first, is about 
16% for the Iwasaki action and about 12% for the plaquette action. Similarly, our final result for the charmed-strange 
meson hyperfine splitting in physical unit is given in Tabl dXlIU where the central value is AMpo\e{D* — Dg), the first 
error is statistical and the second is a systematic error. It is interesting that the second errors for the charmed-strange 
meson hyperfine splitting , which are about 8% for the Iwasaki action and about 7% for the plaquette action, are half 
of those for the charmonium hyperfine splitting. This suggests that the dominant systematic errors come from the 
heavy quarks, so that they are proportional to a number of heavy quarks in the mesons. In Table IXIIII our final result 
for the Dg meson decay constant in physical unit is also presented, where we take fDsi-A4) with Mpdc as the central 
value. The first error is statistical and the second and the third are systematic errors estimated from a difference of 
Id, {^4) between Afpoio and Afkin and a diflference between foA^^) and foA^k) with Mpoic, respectively. Both the 
second and third errors are less than 1% for the plaquette gauge action. For the Iwasaki gauge action, on the other 
hand, the third error is about 5% though less than 1% for the second. Smallness of the third error for the plaquette 
action may be partly due to the use of c^'^"'"'^^. Note that the systematic errors associated with the heavy quark 
action are estimated at one lattice spacing in this paper. Therefore, in future works, it is desirable to study these 
systematic errors by changing the lattice spacing. 

Once the systematic errors are taken into account, our results of the hyperfine splitting for two gauge actions agree 
with each other. For fu^, on the other hand, an agreement is not so excellent: the difference is still 1.5a even if we 
take the systematic error for the Iwasaki action. It could be interesting to see whether the difference diminishes if we 
employ c^"'"'*'^^ for the Iwasaki gauge action. 

D. Comparison with tlie previous results 

In Fig |12l our results of the S-state charmonium hyperfine splitting are compared with a previous result obtained 
by the CP-PACS collaboration using the anisotropic lattice QCD[1^, where the effective speed of light is nonper- 
turbatively adjusted to unity such that M^°^° = AI^"^. Both results are plotted as a function of the lattice spacing 
determined by the Sommcr scale ro = 0.5fm. Our result with the pole mass for the Iwasaki gauge action is consistent 
with the continuum limit of the anisotropic lattice result within the small statistical error, though the kinetic mass 
result is rather large. For the plaquette gauge action, on the other hand, both the pole and kinetic mass results are 
larger than the anisotropic lattice results. The large systematic error due to the pole to kinetic mass difference should 
be eliminated with the use of nonperturbative i' in future calculations. It should be noted that all the results are 
smaller than the experimental value by about 40%. 

Figure [13] shows the comparison of our results of the S-state charmed-strange meson hyperfine splitting with a 
previous result obtained by the UKQCD collaboration using the heavy clover quark action[20]. We observe that all 
the results agree within large statistical errors, though they are smaller than the experimental value by about 10%. 

In FiglU] we compare our results of fo^ with a previous result obtained by the ALPHA collaboration using the 
heavy clover quark and the plaquette gauge actions [2l|. Our results at finite lattice spacing are closer to the ALPHA 
result at the continuum limit than at a similar lattice spacing. This could indicate that fo^ from the RHQ action has 
a good scaling behavior, which should be checked in future scaling studies. We also point out that c^^^^^ for the 
Iwasaki gauge action may reduce the difference between foS^i) and foA^k)- We also leave it to future work. 
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V. CONCLUSION 



We have carried out a first nonperturbative test of the RHQ action focusing on the magnitude of the cutoff errors. 
We investigate the dispersion relation of the pseudo-scalar meson and the space-time symmetry for the pseudo-scalar 
meson decay matrix element. Our results show that the RHQ action has much smaller cutoff errors than the heavy 
clover quark action around the charm quark mass. 

We also investigate the systematic errors due to the cutoff effects for the physical observables. In case of the 
charmonium (charmed-strange) hyperfine splitting, a difference between the results with Af and M*^'" is used to 
estimate the systematic error, which is as large as 16% (8%) for the Iwasaki gauge action and 12% (7%) for the 
plaquette gauge action. For the Dg meson decay constant fo^ i we estimate the systematic error by a difference 
between JdX-^'^) ^^'^ fo^iAk) as well as a difference between Mp°^° and M'^'". The latter is negligible for both gauge 
actions, while the former is about 5% for the Iwasaki gauge action and 0.5% for the plaquette gauge action. 

There are two important subjects for future studies. One is a further improvement of the RHQ action to reduce 
the cutoff effects. In particular, it is rather easy to tune the improvement coefficient i> nonperturbatively, which is 
supposed to eliminate the leading 0{al) errors. This study is under wavp^. The other is the inclusion of light 
dynamical quark effects. It is interesting to investigate whether the deficit in the quenched value for the S-state 
charmonium hyperfine splitting is fully accounted by the sea quark effects. 
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APPENDIX A: RENORMALIZATION FACTORS AND IMPROVEMENT COEFFICIENTS FOR 

MASSIVE QUARKS 

In this appendix we explain how to determine the input parameters for the RHQ action and the axial vector currents 
in our numerical simulation, such as n, the iniprovement coefficients and the renormalization factors together with 
the mean field improvement discussed in Refs.fa, Q. 

The mean field improvement is introduced as the redefinition of link variable U^{x) —>■ uo(C/^(x)/uo) = Uf)U^{x), 
where mq = P^/'' with the averaged plaquette value P in our simulation. The one-loop expression for uq is given by 

"0 = 1- g^^TiviF, (Al) 

where Tmf = 1/8 for the plaquette gauge action and 0.0525664 for the Iwasaki gauge action[23j. 

With the replacement J7^(x) uqU^j_{x) it is natural to introduce the boosted gauge coupling Qq/uq, which is 
related to the MS coupling constant gl^(fi) with the scale /i = 1/a as 

1 cqP + 8ciR n-^Nf ^ ^ ^ 

2 = 0.1006 + 0.03149Ar/ + — ^ log{an) (A2) 



for the Iwasaki gauge action and the 0(a) improved Wilson quark action[2J|, and 

IP 11 - -iVf 

„ , , = ^ - 0.1349 + 0.03149iVf + 1— ^-log(a/i) (A3) 

gjj^if^) 9o 87r2 

for the plaquette gauge action and the 0{a) improved Wilson quark actionflTj. In the following we simply use to 
express 9j-^{fJ. = 1/a)- 

The inverse quark propagator at the leading order without the mean-field improvement is given by 

aS~^ = amo -f g^am^^-* -I- 174 sin(p4a) + i;/ 7/c sin(pfca) -|- ri(l — cos(p4a)) + Ts ^^(1 — cos(pfca)) (A4) 
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where amo is the bare quark mass appeared in the action. Note that we include the one-loop contribution to the 
critical quark mass, g^mi^\ in the leading order. A reason for this will become clear later. The pole mass am^\ 
determined from the zero of the inverse propagator by setting p4 = im'p^ and pk = 0, satisfies 

sinh(am^°^) + rt cosh(ar7i^*'^ ) = am + rt, (A5) 

where am is a shifted quark mass defined by am = amo + ag'^mi^K If we perform the replacement U^{x) — > UQUf^{x) 
in the RHQ action given in Eq.([T]) or Eq.Q, the inverse quark propagator at the leading order with the mean-field 
improvement becomes 

aSq^ — amo + g'^aAm[^^ + 174^0 sin(p4a) + 7fcUo sin(pfca) + rt{l - uo cos(p4a)) + ^(1 - uq cos{pka)) 

k k 

(A6) 

where g'^a/S.m^c^ = g'^am^c^ — (rt + 3rs)(l — uo). Then the pole mass am^^ at the tree-level with the mean-field 
improvement satisfies 

Uo sinh(am^°'') — amo + ag'^Am'^^^ + rt{l — uq cosh(aTO|j"-')) + 3rs{am'^^){l — uo) 

^ am + rtuoil - cosh(am(,°))) + (1 - uo)3{rs{am'j^'^) - 1). (A7) 

Note that the shifted quark mass am is kept equal with and without the mean field improvement. Therefore both mp"^ 
and rfip^'' vanish at am — 0. Since the remaining one-loop correction to the quark mass is multiplicative to m, the pole 
masses in both definitions vanish at am = also at one-loop level. The inclusion of g'^am^c^ or g'^alS.m^c^ at leading 
order is necessary to satisfy this property. Although in this work we follow the mean-field improvement procedure 
given in Sec. 6 of Ref. 8] which does not include the Ato^ correction, the effects on the improvement parameters are 
less than 1%. Eas. (IA5p and (|A7p lead to the following relation that 



foi ~ fol , sinh(aTO^"^) -|- rt(cosh(aTO^°') - 1) + 3(rs(aTO^"-') - 1) 

P = am\,> + {uo- 1) ^( . ,(0). 

cosh(amp ) + rt smii(aTOp ' ) 



amp 



am, 



p 



g'^aAmp, (A8) 



where 



Cprj, sinh(am^"^) + rt(cosh(a7Ti^"^) - 1) + 3(rs(am^"^) - 1) 

aAmp = -— Tmf ,(0) . , , , (o). ' i^^> 

^ cosh[amp ) + rt smh(amp ') 

As a consequence, the quark pole mass is written at the one-loop level as 

amp = am^°^ + g'^am'^p'^ = amf^ + g'^am'^p\ (AlO) 

where afn'p^ = am^p^ + aAmp, and am^^ is the one- loop correction to the pole mass without the mean field 
improvement 0. 

The mean-field improved parameters Zq, Ts, ce and cb are given below with the use of arn^p^ and afn'p^ : 



7(1) ^ 2 

ZqMttiam^p ') - ^g.iatt^o 1+5 + 9 ^Tmf + —^aAmp ) , (All) 



Ka4°)) = + ff^^fi) + g^^%^aAmp, (A12) 



dam^^ 

a (0) 
or: ' 



,,(am(0)) = r^^)+g\m+g^-^aAmp, (A13) 

damp 

1 / (1) q 2 n (0) \ 

c^(amf ) = ^3 ( 1 + - g^'-CT^^ + ^^aAm^ ) , (A14) 

(1) q „2 aJo) 



CB(amW) = 1 + .g2^ - g^^CpT^^ + ^^aAm^ , (A15) 



\ Cr 
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where one- loop corrections, ^^Vatf ^^^^ ^ ^^^^ ^e'' ^^"^ have already been calculated in ref. [1]. We replace 
a perturbative value of uq in the above expressions by uq ~ p^/'^ with P taken from our simulation. We finally 
determine k in terms of m^^ as follows. Using the relation of Eq. (jA7P with rt — 1 



ama + aArric ^ uoe'^"^p - 1 - 3rs(l - uq), (A16) 



the hopping parameter k is given in terms of m'^"' : 

K=- i = - -j^ ^- . (A17) 

2 1 + 3r^ + amo '2 uq^c'^'^p + Sr^) ~ aArric 

With this definition, k becomes Kcrit at the one-loop level for arh^'^ — 0. 

In a similar manner we can derive the renormalization factor and the 0{a) improvement coefficients for the axial- 
vector currents in Eq.®. The matching factor Za,^, from the lattice to the continuum MS scheme is given in Ref.0]: 



^latt 

rz A' /'7(o) / ~(o)n //7(o) / ~(o)n /, : 

" ^If = V ^Q.iatt(«'^pi ) V ^9,iatt(a"ip2 )wo (1-3 

Zj 



A, 



+5 ^tmf + 2 ^(0) . (0) "^'"^i + 2^^;;75r«^"^^^ ' (^^^^ 

^Q.latt '^"'pl '^g.latt '^"'■p2 / 

where A^^ is the one-loop correction to the renormalization factor of ^^[9]. For the 0(a) improvement coefficients, 
on the other hand, we use the expressions of Eq.® suitable for our numerical simulations, which are related to those 
in Ref.[§l as 

(A19) 

afa-^ ) /sinh(aTO(") -|- am-^ ) , (A20) 

(A21) 
(A22) 
(A23) 

where c^^^'^'^^'^""" are calculated as a function of •ifiQ and the superscript PT represents that these parameters are 
defined in Ref.0. Note in particular that a minus sign in the relation (|A22p . A factor 1/ito in Eqs. (|A2ip — (|A23|) is 
due to link variables in the point splitting operators of Eqs. (fTO|) and (fTTj) . In Eq. (|A20p we multiply an extra factor 

[rniP -|- mQ^)/sinh(m^°^ -f mg^) since c^^'^ in Ref.[9|] is a coefficient of [m'q^ -f ■mQ^)q{q)Q{p) while c^^ in Eq.® is 
a coefficient sinh(mq°-' + mQ'')q{q)Q{p). 





2 +,PT 

= 9 , 




2 +.PT, ~ (0 




2 -,PT / 
= 9 /"O, 




2 i.PT / 
= -9 /Uq, 




2 -R.PT / 

= g /uo, 
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TABLE I: Simulation parameters for the RHQ action with the Iwasaki gauge action(upper) and the plaquette gauge ac- 
tion(lower). 





















flavor 


K 


Mps/Mv 


1/ 






— Ctt 1 






0.13295 


0.5567(36) 


1.0 


1.0 


1.50 


1.0 


llglil 


K2 


0.13222 


0.6898(23) 




1 n 




1 n 






0.13138 


0.7734(16) 


1.0 


1.0 


1.50 


1.0 






0.11513 


0.9372(7) 


1 03160 


1.12787 


1 66304 


n 92064 






0.10524 


0.9680(4) 


1 05935 


1.20160 


1.75930 


n ssssq 






0.09455 


0.9813(2) 


1.10040 


1.29777 


1.88490 


0.85628 






0.07841 


0.9901(1) 


1.19159 


1.48857 


2.13426 


0.81003 










plaquett 


e 








flavor 


K 


Mps/Mv 


V 




cs 


UJ — ce/cb 




ni 


0.13449 0.5492(173) 


1.0 


1.0 


1.769 


1.0 


hght 


K2 


0.13373 


0.7088(27) 


1.0 


1.0 


1.769 


1.0 




H3 


0.13298 


0.7837(16) 


1.0 


1.0 


1.769 


1.0 




K4 


0.11456 


0.9345(5) 


1.04161 


1.16034 


2.02423 


0.91709 


heavy 


t5 


0.10190 


0.9728(2) 


1.08301 


1.26201 


2.17790 


0.88345 






0.09495 


0.9808(2) 


1.11284 


1.32840 2.27776 


0.86586 




m 


0.07490 


0.9911(1) 


1.23871 


1.58259 


2.66050 


0.81742 



TABLE II; Simulation parameters for the heavy clover quark action with the Iwasaki gauge action(upper) and the plaquette 
gauge action(lower). 







Iwasaki 








flavor 


K Mps/Mv V rs 


Cb w 


= ce/cb 


hght 


K3 


0.13138 0.7734(16) 1.0 1.0 


1.50 


1.0 




K4 


0.1256 0.9353(5) 1.0 1.0 


1.50 


1.0 


heavy 


ks 


0.1186 0.9717(3) 1.0 1.0 


1.50 


1.0 




K6 


0.1119 0.9836(2) 1.0 1.0 


1.50 


1.0 






plaquette 








flavor 


K Mps/Mv V Ts 




= ce/cb 


light 




0.13298 0.7837(16) 1.0 1.0 1.769 


1.0 




K4 


0.12780 0.9359(8) 1.0 1.0 


1.769 


1.0 


heavy 


K5 


0.11900 0.9770(3) 1.0 1.0 


1.769 


1.0 




He 


0.11480 0.9834(2) 1.0 1.0 


1.769 


1.0 



TABLE III: Smearing parameters. 



action 


Iwasaki+RHQ 
flavor A B 


Iwasaki+CL 
flavor A B 


plaquette+RHQ 
flavor A B 


plaquette+CL 
flavor A B 


light 


Ki 1.28 0.28 
K2 1.25 0.3 
K3 1.25 0.32 


K3 1.25 0.32 


Ki 1.28 0.35 

1.28 0.35 
K3 1.28 0.35 


K3 1.28 0.35 


heavy 


K4 1.25 0.50 
K5 1.25 0.58 
K6 1.25 0.65 
K7 1.25 1.00 


K4 1.25 0.50 
K5 1.25 0.65 
Ke 1.25 0.82 


K4 1.25 0.50 
K5 1.25 0.58 

K6 1.25 0.60 
K7 1.25 0.8 


K4 1.25 0.50 
/t5 1-25 0.65 
K6 1-25 0.65 
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TABLE IV: Fitting range from tmin to tmax for the two-point frinctions. 









Iwasaki 


plaquette 


correlator 


system source t 


mill / ^max 


tmin/ ^max 




H-L 


S 


10/22 


7/17 




H-L 


s 


10/22 


9/22 


(Ppt) 


H-L 


p 


13/22 


14/22 


(Ppt) with IpI ^ 


H-L 


s 


8/20 


9/22 


{VV^} 


H-H 


s 


12/22 


11/21 


(ppt) 


H-H 


s 


12/22 


11/21 


(ppt) 


H-H 


p 


16/22 


18/23 


(PPt) with \p\^0 


H-H 


s 


12/20 


10/22 



TABLE V: Fitting resuhs tor chiral ('xtrapcjlaticni of Mrs f"i' tl"' light-light pseudoscalar meson. 



action 




xVdof 


Iwasaki+RHQ 


0.7583(15) 0.133802(4) 


6.3 


plaquette-|-RHQ 


0.7886(38) 0.135247(9) 


0.7 



T.VBLE \'l: Hopping paraiuotors at the physical pcjiiits. 



action 




Kci J /-([>, Mpole) 


Kc(J/V',-^kin) 




Kc{Ds,Mkin) 


Iwasaki-I-RHQ 
plaquette-|-RHQ 


0.133749(4) 0.132422(4) 
0.135200(9) 0.134026(5) 


0.099414(22) 
0.100593(21) 


0.102362(377) 
0.102610(343) 


0.099640(49) 
0.100669(35) 


0.101528(862) 
0.102402(810) 



TABLE VII: Fitting results for the heavy-heavy hyperfine splitting as a function of the vector meson mass. 



action 


A 


Mpole 

B C xVdof 


A 


Mkin 

B C 


xVdof 


Iwasaki-|-RHQ 
plaquette-l-RHQ 


0.0452(14) 
0.0524(14) 


0.0156(43) 0.0009(26) 0.0018 
-0.0028(33) 0.0081(19) 0.03 


0.0676(23) 
0.0745(20) 


-0.0120(78) 0.0122(50) 
-0.0403(59) 0.0288(34) 


0.001 
0.10 



TABLE VIII: Charmonium and charmed-strange meson hyper-fine splittings in unit of GeV with Mpoie and Mkin- The lattice 
spacing is determined by the Sommer scale ro = 0.5fm. 





Iwasaki plaquette 
Mpole Mkin Mpole Mkin expt. 


AM(J/i/) - r)c) 
AM{D: - Ds) 


0.0728(8) 0.0847(20) 0.0788(6) 0.0877(18) 0.1173(12) 
0.1243(28) 0.1348(65) 0.1261(16) 0.1358(57) 0.1438(4) 



TABLE IX: Fitting results for the heavy-light hyperfine splitting as a function of the pseudoscalar meson mass. 



action 


Mpole 

ABC 


xVdof 


Mk 

A B 


n 

C 


xVdof 


Iwasaki-|-RHQ 
plaquette-|-RHQ 


0.0602(35) 0.0028(17) -0.0628(167) 
0.0545(20) 0.0022(10) -0.0523(121) 


0.13 
0.22 


0.0649(75) 0.0012(41) 
0.0549(70) 0.0049(36) 


-0.0106(337) 
-0.0224(497) 


0.03 
0.61 
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TABLE X: Fitting results for the heavy-light pseudoscalar meson decay constant determined from A4, as a function of the 
pseudoscalar meson mass. 



action 


A 


B 


C D 


xVdof 


Iwasaki+RHQ 
plaquette+RHQ 


0.2204(93) 
0.1285(20) 


-0.1337(113) 
-0.0217(26) 


0.0312(39) 0.1780(199) 
-0.0065(10) 0.1967(140) 


0.43 
2.9 


action 


A 


B 


C D 


xVdof 


Iwasaki+RHQ 
plaquette+RHQ 


0.2275(96) 
0.1326(29) 


-0.1470(127) 
-0.0241(49) 


0.0368(47) 0.1820(210) 
-0.0071(22) 0.1903(131) 


0.33 
1.7 



TABLE XL Fitting results for the heavy-light pseudoscalar meson decay constant determined from Ak as a function of the 
pseudoscalar meson mass. 



action 


A 


A/polo 

BCD 


xVdof 


Iwasaki+RHQ 
plaquette+RHQ 


0.1892(88) 
0.1709(57) 


-0.0955(110) 0.0174(39) 0.1788(190) 
-0.0811(73) 0.0144(27) 0.1338(306) 


0.44 
0.55 


action 


A 


BCD 


xVdof 


Iwasaki+RHQ 
plaquette+RHQ 


0.1975(96) 
0.1746(54) 


-0.1104(129) 0.0235(50) 0.1836(198) 
-0.0868(70) 0.0163(28) 0.1370(311) 


0.29 
0.86 



TABLE XII; Ds meson decay constants in unit of GeV determined from Ak and A4, using iV/poio as well as Mkin. The lattice 
spacing is determined by the Sommer scale ro = 0.5fm. The experimental value for fn^ is take from Ref . [l^ . 





Iwasaki plaquette 

Mpolo Mkin Mpole Mkin expt. 




0.2506(49) 0.2496(48) 0.2291(22) 0.2304(24) 0.282(16)(7) 
0.2373(47) 0.2369(46) 0.2305(31) 0.2304(30) 0.282(16)(7) 



TABLE XIII: The final results of the charmonium hyperfine mass splitting, the charmed-strange meson hyperfine splitting 
and the Ds meson decay constant in unit of GeV. The first error is statistical and the second and the third ones are the cutoff 
errors explained in the text. The lattice spacing is determined by the Sommer scale ro = 0.5fm. 
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FIG. 1: Effective mass plots for the H-L(left) and H-H(right) mesons with zero spatial momentum in case of Kheavy = and 
Kiight = K3 with the Iwasaki gauge action. Circles represent the local source correlators and squares for the smeared source 
correlators. 
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FIG. 8: Same as Fig[7]for the heavy- light axial vector current. 




FIG. 10: Heavy-light meson S-state hyperfine splittings as a function of aPIps with Kught = ks- 
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FIG. 11: Charmed-strange pseudoscalar meson decay constants obtained from the temporal(top) and spatial(bottom) compo- 
nents of the axial vector current as a function of aMps in case of KHght = K3- 
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FIG. 12: Comparison of charmonium S-state hyperfine splitting in physical unit. 
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FIG. 13: Comparison of Ds meson hyperfine splitting in physical unit. The results of Ref.[20| are slightly shifted in a coordinate 
for visibility. 



0.32 



0.3 



0.28 



0.26 



0.24 



0.22 



0.2 



fo [GeV] 



* 


expt. [17] 




• 


RHQ+plq, 




o 


RHQ+plq, 




■ 


RHQ+Iwa, A^ 




□ 


RHQ+Iwa, A^ 




A 


CL+plq, A4 [20] 





i 



i 



_L 



0.02 0.04 0.06 

a[fm] 



0.08 



0.1 



0.12 



FIG. 14: Comparison of Ds meson decay constant in physical unit. 



